Wetting to Non-wetting Transition in Sodium-Coated Ceo 
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Based on ab initio and density-functional theory calculations, an empirical potential is proposed 
to model the interaction between a fullerene molecule and many sodium atoms. This model predicts 
homogeneous coverage of Ceo below 8 Na atoms, and a progressive droplet formation above this 
size. The effects of ionization, temperature, and external electric field indicate that the various, and 
apparently contradictory, experimental results can indeed be put into agreement. 
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Fullerene molecules are now commonly used as build- 
ing blocks of complex materials having unusual physical 
and chemical properties H. The possibilities offered by 
molecules such as C60 or C70 in terms of electric or op- 
tical devices have laid ground for the rapidly expanding 
nanotechnology area. The interaction between fullerenes 
and alkali atoms has attracted a lot of attention after 
the discovery of superconductivity in the K 3 C 60 solid 
Sodium-Cgo compounds have also received much atten- 
tion H H . More recently, several groups investigated ex- 
perimentally the gas phase properties of exohedral alkali- 
C60 compounds || g, f?J , to further predict the solid state 
of alkali-doped fullerite. 

Using mass spectrometry, Martin and coworkers || 
inferred from the most stable "magic numbers" that 
sodium covers the Ceo molecule in a continuous and 
homogeneous fashion, and that metallic bonding starts 
above 6 metal atoms. Palpant and coworkers de- 
duced from photoelectron spectroscopy measurements 
that coating proceeds by trimers units rather than single 
atoms. They also estimated that metallic bonding only 
appears at n = 13 sodium atoms. Very recently, Dugourd 
and coworkers measured the electric polarizability of 
Na„Cgo clusters in the range 1 < n < 34. They con- 
cluded that a sodium droplet is formed on the surface of 
the fullerene. These apparent contradictory results are 
surprising, as far as they rely on similar experimental 
conditions but different interpretations. 

The situation is in fact even more intricate due to the 
variety of theoretical conclusions on the very same sys- 
tems H fuj [ll], 0. With the exception of the study by 
Rubio and coworkers which assumes complete wet- 
ting of Ceo by sodium in a continuous two-shell jellium 
description, the electronic structure calculations by Hira 
and Ray [fl], and by Hamamoto et al. Jl2| reach 
different conclusions as to whether the Na2 molecule re- 
mains as a dimer or dissociates into atoms, which locate 
on opposite sides of Cgo- At the unrestricted Hartree- 
Fock level, Hira and Ray find that charge transfer is neg- 
ligible in NaCgo, in apparent contradiction with experi- 
ments 0, D , and also that Na2 remains as a dimer loosely 
bound to Cgo- On the other hand, the more realistic 
density-functional theory (DFT) calculations performed 



by Hamamoto and coworkers |12| at the local density 
approximation level tend to favor regular coating, pref- 
erentially by trimers, in agreement with the picture of 
Palpant et al. (7). However, all these ab initio investi- 
gations could only focus on very limited sets of chosen 
geometries, simply because global optimization on such 
complex energy landscapes is far beyond the current pos- 
sibilities of fully quantum mechanical approaches. In ad- 
dition, statistical and temperature effects still remain a 
major limitation in large scale simulations based on first 
principles. 

To overcome the above difficulties, we have constructed 
an empirical, ab initio based atomistic model allowing 
systematic static and dynamical investigations in a large 
size range. Briefly, given the geometry R = {r^r^} of 
the Na n Cgo^ system, with € Na„ and r'j 6 Ceo, the 
total potential energy is written as 



V(R) = TWl-z) + ^CeoK') + Mntcr(R) + Vboul(R), 



(1) 



where VNa„, Vc 60 , Vinter, and Vboui are the energies of 
sodium, Ceo, covalent interactions between Na and Cgo, 
and Coulomb interactions, respectively. Metallic forces in 
sodium are described by a simple many-body potential in 
the second-moment approximation of the electronic den- 
sity of states in the tight-binding theory [|l3| . Covalent 
forces in Cqo were most often neglected, as the temper- 
atures of interest are low enough (< 300 K) to keep the 
fullerene rigid. We occasionaly checked the rigidity of Cgo 
in some optimal structures with the Tersoff potential [ 14j . 
The covalent interaction between sodium and carbon 
atoms was modelled as a pairwise repulsion term taken 
from ab initio calculations of the 4 n state of the NaC 
molecule, and expressed as ^N°£ alont (r) = D exp(-ar). 
To account for charge transfer, and to quantify the ex- 
tent of ionic interaction in the system, fluctuating charges 
Q = {lii with qi G Na„ and q'j € Cgo, were added to 
all atoms in the system. For any geometry R, the charges 
are determined in order to minimize the global Coulomb 
interaction Vcoui , which is equivalent to equalizing the ef- 
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fective electronegativities in the sense of Sanderson |15|] : 

Vq(R) = E £Na * + ^ Na ^ + E ec ^ + ^c9; 2 

i j 

+ E J H'9iQi' + E ■'■J'' 1 ',' 1 '. 

*<*' i<i' 

+ E -m<'/' + A ( E « + E ■ - q ) -(2) 



In this equation, i and j label sodium and carbon atoms, 
respectively. £Na and £q are the respective electronega- 
tivities of sodium and carbon, i?Na and Hq their 'hard- 
nesses', Jij the Coulomb interactions. The hardnesses 
parameters correspond to the on-site repulsion in the 
Hubbard model. The latter quantities were extracted 
from ab initio calculations of the Coulomb integrals. Due 



to the finite derealization of the electrons, does not 
diverge at short distances. We have used the simple ex- 

1/2 

pression Jy(r) = [r 2 + H~ 2 exp(— jijr 2 )] to model 
the Coulomb parts of the Na-Na, C-C, and Na-C inter- 
actions. The hardnesses parameters H^ a and He are the 
r — > limits of JNa-Na and Jc-G- Finally, the Lagrange 
multiplier A in (||) ensures that the system carries a global 
charge Q. This electrostatic model was further improved 
by including 90 fixed charges, 60 on each carbon site 
with value Sq, 30 on the middle of each C=C bond with 
value —2Sq. This was previously used by Schelkacheva 
Tareyeva in a study of bulk Cgo fl6| ■ 

The dynamical use of the full potential energy func- 
tion is greatly facilitated by considering an extended La- 
grangian where the fluctuating charges are treated as in- 
dependent variabl es |17| , in a way similar to the Car- 
Parrinello scheme |18[ . This empirical model has 15 pa- 
rameters, namely the 5 parameters of the SMA potential, 
the 2 parameters of the covalent Na-C interaction, the 
electronegativity difference ec — £Na, the 6 parameters 
of the Coulomb interactions (including hardnesses), and 
the fixed charge Sq. These parameters were obtained 
by minimizing an error function \ 2 to reproduce sev- 
eral properties independently evaluated by ab initio and 
DFT calculations [gj| with the B3LYP nonlocal func- 
tional ^0|. These properties are the charge transferred 
from sodium (0.87e) and the electric dipole (14.5 D) in 
NaCeo, the energy difference between the two Na2C6o iso- 
mers with sodium atoms on adjacent or opposite hexag- 
onal sites (AE = —0.35 eV), and the binding energy 
and equilibrium distance in Na2 (resp. D e =0.73 eV and 
R e =3.08 A). It should be noticed that these theoretical 
values agree well with the experimental data available 
H . We also added to the error function a penalty term 
to reproduce the experimental electric polarizability of 
the C 60 molecule (76.5 A 3 |l|). 

The lowest energy structures of Na n C6o clusters were 
determined using a variant of the basin-hopping, or 
Monte Carlo + minimization algorithm p2| , in which 
each sodium atom is offered a probability to rotate freely 



over the Cgo surface, in addition to the usual random 
displacements moves. Some of the optimal structures are 
represented in Fig. |l| for the sizes n — 4, 7, 8, 12, and 20. 
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FIG. 1: Some optimal structures of Na n C60 clusters, (a) 
n = 4; (b) n = 7; (c) n = 8; (d) n = 12; (e) n = 20. The 
lowest energy structure for Na2C60 within an electric field 
E = 2.55 x 10 s V.m" 1 if represented in (f). 

Up to seven sodium atoms, the empirical model predicts 
"wetted" structures, where each alkali atom stays over 
an hexagonal site. In this regime, the Coulomb repulsion 
between the positive charges carried by the sodium atoms 
is minimized by keeping the atoms as far away as possible 
from each other. At the stage of 7 atoms, there is only 
little space left on the Cgo surface to add an extra eighth 
isolated atom, and the energy gained in creating Na-Na 
bonds gets larger than the electrostatic penalty. This em- 
bryonic droplet further grows as new sodium atoms are 
incorporated, and progressively spreads over the Ceo sur- 
face. Growth occasionally proceeds by removing one of 
the remaining isolated atoms on the side opposite to the 
droplet, and may eventually end in a single, big droplet. 
The largest system considered here, Na3oCgo, however, 
still shows one main droplet and two isolated atoms. 

To provide a more complete interpretation of the ex- 
periments, we have considered the separate effects of tem- 
perature, ionization, and a possible external electric field 
on the transition between coated and segregated mor- 
phologies. Equilibrium dynamics at constant tempera- 
ture was achieved by supplementing the extended La- 
grangian with two sets of Nose- Hoover chain thermostats 
|23[ , one set at the desired ionic temperature T, one set 
to keep the average kinetic energy of the charges to a low 
value T* = T/100. Thermostating the charges also pre- 
vents the fictitious dynamics to diverge too much from 
the adiabatic, Born-Oppenheimer type dynamics [^4| . 
Diffusion constants were estimated from the slopes of the 
average mean square atomic displacements versus time in 
a set of 500 simulations at various temperatures in the 
range 25 K< T <300 K, with increments 25 K. In Fig. |, 
these diffusion constants are represented for the sizes 4, 
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8, 12, and 20 as a function of 1/T in an Arrhenius plot. 
These plots provide us with estimates of the correspond- 
ing activation barriers, for which we find A ~ 400 K for 
all four sizes. At room temperature, sodium atoms thus 
show a significant mobility over the Cgo surface. This 
agrees with the recent experimental findings on KCgo 
p5[ . To quantify the effect of temperature on the mor- 
phology in Na„C6o, we have calculated the average size 
of the largest sodium fragment. A fragment is here de- 
fined as a set of connected atoms, i.e. having at least 
one neighbor at a distance shorter than 8 bohr. The 
same geometrical observable was also computed at T = 
for neutral or singly charged molecules. In the case of 
charged systems, global optimization of Na„C6o + and 
Na n C6o~ without electric field was performed for sizes 
near the crossover, n < 12. Lastly, we have compared 
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FIG. 2: Arrhenius plots of the diffusion constants of some 
Na n Ceo clusters, for n = 4, 8, 12, 20 (no field), and 6 within 
an electric field. 

the results without electric field and those with a homo- 
geneous field of magnitude E — 2.55 x 10 8 V.m -1 relevant 
to the experiments of Dugourd et al. § , and at zero tem- 
perature. The corresponding data were averaged on 30 
different orientations of the field, relative to the Cgo ge- 
ometry. Global optimization was carried out for each of 
these orientations. The variations of the largest fragment 
size (N) against the total number of sodium atoms n are 
plotted in Fig. |. The rise of (N) at the crossover size n 
marks the onset of droplet formation. As temperature 
increases, isolated atoms are more likely to meet either 
with each other or with the preexisting droplet. This 
lowers n* by a visible amount. However, 300 K is a rela- 
tively high temperature for sodium clusters lying usually 
above their melting points p3] . Thus, external bonds are 
often broken in the course of the dynamics, which results 
in lowering (N) at larger n. Ionizing the cluster can also 
strongly affect the crossover size n*: negatively charged 
molecules are more stable in segregated form due to a 
much lower Coulomb repulsion of the alkali atoms, and 
positively charged molecules are less stable and remain 
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FIG. 3: Size of the largest fragment of sodium atoms in 
Na^Ceo clusters in various conditions, and for neutral or 
charged molecules. 

At the experimental value, a finite electric field dis- 
favors regular location of the charges over the nearly 
spherical Cqq surface. This yields more frequent metallic 
bonds. In addition, the lowest energy structures found 
for all field directions can also change significantly. For 
instance, as seen in Fig. |l], two sodium atoms get closer 
to one another, yet they still do not form a bond at K. 

Electric susceptibilities x have been estimated with the 
present model using the high-temperature approximation 
to the Langevin theory for rigid dipoles undergoing ori- 
cntational thermal fluctuations: x — a + H 2 /ik^T . In 
this formula a is the static polarizability, fee is the Boltz- 
mann constant, and fi the electric dipole. The variations 
of x with cluster size n are shown in Fig. [i[ They qual- 




FIG. 4: Electric susceptibilities at 300 K. Inset: electric dipole 
moment and comparison with DFT calculations. 

itatively resemble the experimental data of Dugourd et 
al. Ipl , but lie somewhat lower than the actual measured 
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values. This may first come from the parameterization. 
We tested this hypothesis by comparing the predictions 
of the present empirical model to single-point DFT cal- 
culations performed at the same geometry. The electric 
dipoles given in the inset of Fig. show good agreement 
between the empirical and electronic structure estimates. 
We also note that the strong variations in the dipole re- 
flect the changes in the molecular symmetry. It is also 
worth pointing out that the model predicts an energet- 
ical ordering of the isomers essentially similar to that 
found by Hamamoto et al. 12 for the very restricted 



set of structures considered by these authors. Therefore, 
the quantitative difference between our results and the 
experimental data more likely comes from the Langevin 
formula, strictly valid for rigid dipoles. This approxima- 
tion may be insufficient to describe the dipole fluctua- 
tions associated with the floppy character of the metal 
droplet. However, regarding the difference between the 
present results and the data obtained in Ref. ||, adding 
several isolated atoms on the opposite side of the growing 
droplet decreases the dipole moment, hence the suscep- 
tibility. A fully dynamical simulation would probably be 
required to get a more quantitative agreement. 

The present investigation shows that the appearance 
of metallic bonding corresponds to seeding a droplet, 
which progressively grows and captures remaining iso- 
lated atoms. This picture partly reconciles the appar- 
ently contradictory experimental interpretations |T . . 
The empirical potential was fitted to reproduce electronic 
structure calculations as well as some experimental data. 
By allowing large scale simulations and structural sam- 
pling unaccessible to quantum mechanical studies, we 
could estimate the crossover size between wetted and 
nonwetted morphologies to be located near 8 atoms. This 
is also the size range where the metallic transition can be 
estimated to occur. We have also been able to investigate 



how the crossover size depends on various effects such as 
charge or temperature. 

The current limitations of the present empirical model 
are mostly due to the lack of an explicit account of elec- 
tronic structure. Even though Coulomb repulsion should 
play an important role when only few sodium atoms are 
effectively charged but close to each other, shell closing 
can enhance stability in a local way. In particular, this 
could help in explaining why charged trimers remain sta- 
ble. Such quantum effects are obviously beyond the fluc- 
tuating charges model, especially with its current train- 
ing on purely neutral molecules. A next step could be a 
simple quantum tight-binding approach, which would in- 
clude both covalent and charge transfer effects. However, 
a realistic tight- binding Hamiltonian for the Na n Cgo sys- 
tem would require 240 + n electronic states (4 s + p elec- 
trons for each carbon atom, 1 s valence for each sodium 
atom). This is still too heavy for exhaustive sampling 
and large scale simulations. In addition, precise quan- 
tum effects are essentially expected at small sizes, and 
the role of a single charge on the structure and stability 
should decrease as the sodium droplet grows. 

For large sizes, it would be interesting to study in more 
details how the physical and chemical properties of the 
sodium droplet are modified by the interaction with the 
Ceo molecule. The influence of the Cqq deformation and 
vibration on the metal cluster dynamics could also be 
an important issue. Use of the present model to treat 
other metals could help in understanding the differences 
observed in the experiments j?], [|, ^7). Further exten- 
sions, including the treatment of endohedral or exohedral 
fullerenes, nanotubes or surfaces with materials possibly 
other than metals [^8], ^9|, could also be carried out. 

The density functional calculations were performed us- 
ing the GAUSSIAN98 software package @. We thank 
CALMIP for a generous allocation of computer resources. 
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